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ABSTRACT 


Suppose  that  two  discrete  time  series  {x^} 


and  {y  } satisfy  the 


model 


y = V h X + 
7t  s s t-s 


Z. 


ax  ,x 
s s t-s-d  t-s 


+ c 


t* 


where  {x^}  and  {e^}  are  mutually  independent  Gaussian  processes.  The 

lag  d,  the  weights  {at},  and  the  spectrum  of  the  error  process  are 

unknown.  If  ) a = 0 or  c (d)  = 0,  where  c (t)  denotes  the  covariance 

function  of  ix^},  then  the  nonlinear  effect  is  not  detectable  in  the 

spectrum  S (oj)  of  {y  },  or  in  the  cross  spectrum  S (w) . The  cross 
y t xy 

bispectrum,  denoted  (o>^  .u^) , ^as  a "bidden  periodicity"  as  a 

function  of  x = + u^.  The  period  is  proportional  to  the  lag  d. 

A simple  smoothed  cross  bispectrum  estimator  is  defined  in  the 

paper.  This  estimator,  denoted  B(tc),  is  a function  of  the  discrete 

Fourier  transforms  of  the  samples  of  the  {x^}  and  { y ^ } processes. 

Assume  that  S (co)  is  constant  in  a band  (0,2uf  ).  Let  d denote  the 
x o 

value  of  6 in  the  range  0 < 6 < T/2  which  maximizes  the  periodogram 

where  T is  the  largest  integer  such  that  T /T  < f and  x = 4xk/T. 
o o — o k 

The  periodogram  has  a peak  of  order  °(Tq)  at  6=d  against  a background 
or  order  0(1).  A heuristic  approach  is  used  to  develop  the  asymptotic 
properties  of  d.  The  estimation  procedure  is  extended  to  a model 


with  multiple  lags. 


Estimating  the  Lag  in  a Nonlinear  Time  Series  Model 
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I n tr od uct Ion 

Suppose  that  two  discrete-time  real  stochastic  processes  {x^}  a .d  {yt) 
satisfy  the  time  invariant  linear  model 

OD 

y =L  h x,  + tfi  (i) 

where  { } is  a zero  mean  noise  process  that  is  assumed  to  be  independent 
of  {xt>.  In  time  series  terminology  (Brillinger  [1]),  lEyt)  is  the 
output  of  a linear  filter  whose  impulse  response  is  { h ^ } , and  whose 
transfer  function  is  the  Fourier  transform 

»(<•>>  =£t=_coV'la)t 


The  filter  is  called  realizable  if  h ^ = h_^  = ...  = 0.  In  economics, 
the  model  with  a realizable  filter  is  called  a distributed  lag  relation- 
ship between  the  two  time  series.  The  principal  domain  for  H(w)  is  0<u ><n. 

Assume  that^Jh  | < °°,  and  that  {xfc}  and  { e } are  real  stationary 
Gaussian  processes.  Then  {y  } is  stationary  and  Gaussian.  Suppose, 
however,  that  the  impulse  response  is  a function  of  the  xf  with  a lag. 

In  this  case  the  process  { y ^ } is  no  longer  Gaussian.  To  be  specific 
about  the  relationship  between  the  filter  weights  and  the  x^  values, 
suppose  that  hg  in  (1)  is  replaced  by 


hs({xt}) 


h + a x , 
s s t-s-d 


(2) 


where  {a^}  are  unknown  weights  and  the  delay  d is  an  unknown  parameter. 
Such  a model  approximates  the  response  of  an  amplifier  whose  linearity 
is  destroyed  by  feedback  inductances  in  its  circuit  caused  by  the  input 
signal. 
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When  h.({x  })  replaces  h^  in  (1),  we  have 

y = 2_.  h xt  + ).  a x.  . x + c . 
} t s s t-s  ^*s  s t-s-d  t-s  t 


The  term  x , x is  the  t-sUl  value  of  the  d+1  st  lag  process  of 
t-s-d  t-s 

{ xfc } (Nelson  and  Van  Ness  [5]).  If  the  lag  d were  known,  the  investi- 
gator can  estima te  as  by  least  squares.  When  the  lag  is  unknown,  the 
estimation  problem  is  non-trivial.  Kedem  [4]  presents  an  estimator  of 
the  lag  using  a measure  called  lagged  coherence,  but  even  the  asymptotic 
properties  of  his  estimator  are  unknown.  Regardless  of  the  procedure 
which  is  used  to  estimate  d,  the  analyst  who  is  aware  that  the  parameters 
of  his  impulse  response  are  functions  of  the  x^'s  is  better  off  than 
one  who  analyses  the  data  or  makes  predictions  based  on  a misspecif ied 
stationary  model  (1).  This  paper  presents  a method  for  estimating  d 
which  is  easy  to  compute  from  the  discrete  Fourier  transform  of  a sample 
of  the  two  series.  The  estimator  uses  statistics  which  a frequency 
domain  oriented  analyst  would  use  in  the  course  of  building  a linear 
predictor  of  y^.  The  estimate  of  d is  a diagnostic  tool  for  an  analyst 
who  is  skeptical  about  the  linearity  of  his  model. 
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XfV  a = 0 or  c (d)  = 0,  then  the  analyst  gets  no  information  about  the 
s s x 

existence  of  the  nonlinear  component  from  the  sample  mean  of  {y^}.  For 

the  rest  of  this  paper  assume  that  c (d)  a =0.  Also  assume  that 

E | T |c2(t  ) < ”. 
i 1 ' x 

Since  {x  } is  Gaussian,  Ex  x x =0  for  all  t , and  t~. 

t 1 1 t_2  3 A ^ j 

Thus  from  (3)  it  is  easy  to  show  that  the  spectrum  of  {y^}  is 


S^(co)  = | H(w)  | 2Sx(w)  + | A(o>)  | 2Sz(u)  + (m)  , 


where  A(w)  = 


a e 
s s 


-lios 


0 < u)  < u 


S (w) 
x 


c (t)( 


-iwT 


(5) 


(6) 


is  the  spectrum  of  {x^},  is  the  spectrum  of  the  error  process, 

and  S (u>)  is  the  spectrum  of  the  d+1  st  lag  process  Z .(d)  = x ,x  . 
z t t“U  L 

By  the  Gaussian  nature  of  {xtJ,  it  follows  that 


S (oj)  (c2(t)  + c (t-d)c  (t+d)]e 


-1U)T 


T X 


(7) 


Suppose  that  the  analyst  knows  the  transfer  function  H(w)  but  not 

the  shape  of  S^Cw).  He  then  cannot  detect  the  presence  of  the  nonlinear 

term  by  estimating  the  spectrum  of  {y^},  or  the  spectrum  of  the  residuals 

process  {y^  - ?*l^h^x^  namely  |A(w)  | 2Sz(u>)  + S^(w).  On  the  other  hand, 

suppose  that  the  analyst  knows  that  is  a constant  for  all  m,  i.e. 

{e  } is  white.  Then  he  can  identify  the  presense  of  the  nonlinearity 

from  the  residuals  spectrum  provided  that  | A (uj)  | is  not  flat. 

This  will  not  be  true,  however,  if  {xf)  is  white  and  A is  flat  since  it 

can  be  seen  from  (7)  that  S is  flat  when  {x  } is  white.  In  other  words, 

z t 


an  analyst  who  uses  a white  input  to  test  a flat  filter  cannot  detect 


4 


nonlinearity  in  the  filter  from  I lie  spectrum  of  the  output  even  if  he 
knows  that  the  error  process  is  white.  There  is  usually  no  compelling 
reason,  however,  to  assume  that  the  noise  is  white  even  if  the  filter 
is  linear. 

Moving  on  to  the  cross  spectrum,  it  follows  from  (3)  that 

Ex  y = 1 h c (t-s),  and  thus  the  cross  spectrum,  denoted  S (to) , is 
t-r7  tssx  K xy 


/ 

'I 

I 
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Thus  from  (8)  and  (9),  the  cross  bispectrum  has  the  form 


S 

xy 


(u>  jWj)  = ACid^ 


+ w„)S  (u>  )S  (d)  ) (e  iu)ld 
2x1x2 


+ 


-iw 

e 


2 


d 


), 


where 


A(w) 


-ldlS 

a e 
1 s 


(ID 


(12) 


The  principal  domain  for  this  cross  bispectrum  is  {0<Wj<_u,  - to 

and  ) . Rather  than  working  in  this  domain,  consider  the  transfor- 

mation of  variables 

tc  — + q>2  » n = (13) 

The  principal  domain  of  the  cross  bispectrum  as  a function  of  this 
variables  is  the  triangular  set  D = {0<k<2tt,  CKti<_2ti:  k + The 

cross  bispectrum  is 

S (x,n)  = S (uOS  (w»)A(ic)cos  (nd/2)e  l(K<*/2)^  (14) 

xy  x i.  x z. 

K—  n 

where  = — — and  ^ • Equation  (14)  will  be  the  cornerstone 

for  the  analysis  of  the  estimator  of  d which  is  given  in  the  next 
section. 


2.  Estimating  the  Delay 

Given  a sample  * - 0 T-l),  compute  the  discrete  Fourier 

transforms 


dx(u>) 

T-l  -iwt 

- Et=o  V 

(15) 

dy(u>) 

- IT'i  y e'iwt 
t=0  7t 

(16) 
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, , . 2i  4i  fin  , A 7i  k , 4 1 

for  the  discrete  frequencies  — , — , — , . . it.  Let  = -y-  and  nj  = — 

denote  the  discrete  frequency  values  for  k = + u>  and  q = cj^  - ■ 

From  expression  (4.3.15)  in  Brillinger  [1],  it  follows  that 

ET“1dx(a)1)dx(W2)d*(oJl  + u2)  = S (ulfu2)  + 0(T_1 ) (17) 

where  the  star  denotes  taking  the  complex  conjugate,  i.e.  d*(w)  = d^(-(u). 
Also  from  expression  (4.3.15),  or  from  Rosenblatt  and  Van  Ness  [6], 

Var  T 1dx(u)1)dx(a)2)d*(to1  + w2)  = TS^o^)  Sx(w2)  Sy  (co^  + cy)  + 0(T  1)  (18) 

in  the  interior  of  the  principal  domain.  Moreover,  d (oi  )d  (iii„)d*(w  + w„) 

x l x z y 2 

and  dx(w^)dx(w^)d*((jj^  + w^)  are  asymptotically  independent  as 
T -*■  “>  if  ujj  u>2  or  a>2  t 

In  order  to  simplify  the  asymptotics,  assume  that  Sx  is  a constant 

for  w in  the  band  0<u)<2nf  for  some  f . Since  the  spectrum  can  be 

o o 

consistently  estimated  from  the  sample  as  T + “ using  one  of  several 
related  estimators  (see  the  review  by  Hinich  and  Clay  [2]),  the  value  of 
the  spectrum  in  this  band  can  be  assumed  to  be  known  and  equal  to  one. 

Let  f^  be  a frequency  such  that  f ^ < 1 - f^,  and  define  to  be 

the  largest  integer  such  that  T^/T  <_  f^.  As  T + “,  T^/T  f^.  Consider 
the  smoothed  cross  bispectrum  estimator  defined  by  the  expression 
T 

B(k)  « T"2lji1dx[(lc+nj)/2]dx[((c-nj)/2]d*(<).  (19) 

The  sum  in  (19)  keeps  rij  in  the  principal  domain  for  each  K<2iifo,  since 
ic+Hj  <2i»f  0 + 2irT^/T<2ir.  Since 


-1  ^"l  ft  1 -1 

T Zj_1co8(njd/2)  - 1 cos (2nxd) dx  + 0(T  ), 

* r\ 


it  follows  from  (14),  (17),  and  (18)  that 


7 


H(k  ) = (2nd)  2A(  sin  2nf  1 ( ,ckd/2)  + i^,  (20) 

By  the  central  limit  theorem,  the  error  term  u^  is  approximately  complex 

normal  N(0,f  S (k  ))  for  large  T.  Moreover,  u and  u , are  approxi ma tely 
1 y k k k 

independent  if  k f k'. 

Equation  (20)  shows  that  Bix^)  contains  a "hidden  periodicity" 
in  k with  period  T/d.  There  is  no  aliasing  as  long  as  d < T/2.  This 
observation  suggests  the  following  estimator  of  the  lag  d.  Let  d denote 
the  value  of  6 in  the  range  0 < 6 < T/2  which  maximizes  the  periodogram 

T_1K°B(k  )e1(Kk6/2)  |2,  (21) 

O K 1 K. 

where  T is  the  largest  integer  such  that  T /T  < f . As  long  as 
o o o 

f^  j d/2,  the  periodogram  (21)  has  a peak  of  order  0(Tq)  at  6=d  against 
a background  of  order  0(1),  and  thus  d is  a consistent  estimator  of  d. 

In  order  to  obtain  the  asymptotic  variance  of  d,  assume  that  S^(ic)  is  constant 
in  the  band  0<K<2tit Thus  the  errors  u^  in  expression  (20)  are 
homoskedastic . 


If  the  u^  variates  were  independent,  the  asymptotic  properties  of 

d are  given  by  Walker  (7].  His  proof  holds  for  this  model  where  the 

3/2  " 

u^'s  are  asymptotically  independent  as  T ■*  <*>.  Thus  Tq  (d-d)  is 

2 ~2 

asymptotically  normal  wi,th  mean  zero  and  variance  6f^(2nd)  S^(|A|sin  2uf^d) 

Since  T /T*f  as  T-*“,  it  then  follows  that  for  large  T, 
o o 

6f  1 (2?id)  2S 

Mean  Squared  Error  of  d *=  — r - ^ . (22) 

T f |a|  sin  2rf,d 
o 1 1 1 

Note  that  if  f^  is  near  zero,  the  asymptotic  variance  of  d is  0(f^  ). 

Several  different  f^  values  should  be  used  in  practise  in  order  to  avoid  the 


possibility  of  setting  f^  near  d/2. 


8 


This  theoret ical  analysis  shows  that 
-lA"lo  . . i(*  6/2)|2 

1 H=lB(Ko)e  k 1 

for  k in  the  region  where  A(k. ) is  slowly  varying  provides  a simple 
and  graphical  method  for  testing  for  the  presence  of  the  nonlinear 
lagged  term  in  the  model  relating  {y  } or  {x^}.  To  give  a name  to  this 
procedure,  call  this  periodogram  the  Lagstrum. 


3 . Mul t iple  Delays 
Suppose  that 


y,  - £ _ h x + L . £ a „ x x + e , 

t s-o  s t-s  Z- 1 s=o  Zs  t-s-d„  t-s  t 


(23) 


where  d ,...,d  are  unknown  delays.  Assume  that  £ a = 0 for  each  Z,  or 
L S s 

c (dp)  = 0 for  Z = 1,...,L.  From  the  previous  analysis  it  follows  that 

X -C 


the 


cross  bispectrum  for  < = + w^,n  = is 

sxy(|C»r>)  = Sx(w1)Sx(w2)Z^,=  1A^(ic)cos(nd£/2)e  (24) 


Assume  that  S and  each  A,  are  flat  in  the  band  (0,f  ).  Then  the 
x Z o 

large  sample  Lagstrum  will  have  L peaks  of  the  order  0(Tq)  at 
6 = d^,(£  = 1,...,L)  against  a background  of  0(1).  There  will  be  no 
aliasing  as  long  as  T >2d^  for  each  Z. 

Once  again  suppose  that  is  also  flat  in  the  band.  As  is  shown 
by  Hinich  and  Shaman  [3],  the  d ,...,d  which  maximize  the  Lagstrum 

1 Ij 

are  maximum  likelihood  estimators  of  the  d^,  ignoring  the  residual 
correlation  between  the  error  terms.  The  Lagstrum  analysis  extends 
to  the  multiple  lag  model  as  long  as  T is  large  and  A^(w)  is  slowly 
varying  in  the  band  (O.f^)  for  each  Z * 1,...,L. 
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